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We investigate the tiiermal entanglement of indirectly interacting two spins through other spins, 
that is, two spins at the ends of a spin chain. We maximize it by tuning the local fields on the two 
spins to obtain the maximized entanglement. We prove that if the two spins are separated by two 
sites or more, there is a critical temperature above which the maximized entanglement vanishes. We 
numerically calculate the maximized entanglement in three-spin chains and four-spin chains. We 
discover that the maximizing local fields on the spins 1 and 2 have asymmetric forms, which implies 
that the asymmetry of the two spins essentially contributes to the entanglement enhancement. 
In the three-spin chains, we explain this enhancement due to the asymmetry qualitatively and 
\ quantitatively in terms of the magnons. In XX and XY four-spin chains, we find that the critical 

temperature shows qualitatively different behavior depending on the conservation of the angular 
' momentum in the z direction. 

(N 

I. INTRODUCTION 

' The quantum entanglement [1, 2] is an essential resource for quantum information processing [3, 4] and hence it 
, is very important to generate the entanglement efficiently. The entanglement in thermal equilibrium (the thermal 
entanglement [5]) has been researched energetically [6-10] because the destruction of the entanglement by thermal 
] fluctuation extremely restricts the entanglement generation. It is therefore of great importance to know how much 
entanglement can be generated in thermal equilibrium [11-13]. We thereby consider here the maximum value of the 
entanglement in thermal equilibrium. In Rcf. [14], it was shown for two spins interacting directly with each other, how 
Ch . much entanglement we can enhance between the spins by tuning the local fields on the spins. In the present paper, 
we will extend the maximization problem to multipartite spin systems; we study the enhancement of the bipartite 
entanglement due to the local fields in the system where the two spins indirectly interact with each other via other 
spins. The main research targets are the following issues: 
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CS) ' • Is it possible to generate the entanglement at high temperatures in the indirectly interacting spins? 

■ • What are the main factors which make it possible for the local fields to enhance the entanglement? 

■ 

' For the directly interacting two spins [14], the answers for these questions have been obtained: the maximized 
fS| . entanglement decays as l/(nnr) as the temperature increases; the entanglement enhancement by the local fields 
' is brought about by the suppression of the thermal fiuctuation, or in other words, the increase of the purity. Our 
, purpose in the present paper is to answer these questions qualitatively and quantitatively for indirectly interacting 
psj ' spins. 

The possibility of the enhancement of the entanglement with the external fields has been shown in many papers [6, 
11-13, 15-17]. However, there have been few reports about the entanglement maximization of specific two spins. 
Moreover, little is known about the mechanism of the entanglement enhancement due to the external fields or the 
condition for the indirect interaction to generate the entanglement. By making these problems clear, we will be able 
to know the fundamental limits of the entanglement generation in thermal equilibrium. 

In our research, we define the maximized entanglement as the maximum value of the entanglement between the 
particular two spins under the condition that we can arbitrarily tune the local fields only on these two spins. We 
regard the other parameters such as interaction parameters as fixed. If the maximized entanglement is equal to zero, 
the entanglement generation is impossible for any values of the local fields. The forms of the maximizing local fields 
reflect basic properties of the entanglement enhancement, and thereby help us solve the above two problems. 

We study general properties of the maximized entanglement between the focused two spins which are connected by 
a spin chain. Our main results arc the following two: 

• At high temperatures, we prove that the maximized entanglement is always equal to zero between the two 
spins which are separated by two or more spins. In other words, above a critical temperature Tc, we can never 
generate the entanglement between spins far apart for any values of the local fields. It is well known [5] that the 
entanglement vanishes above a certain temperature because of thermal disturbance, but if we can control the 
local fields on the focused spins, the increase of the purity owing to the local fields may suppress the effect of 
thermal fluctuations. Indeed, in two and three spin chains, the non-zero entanglement can be achieved at any 
temperatures by increasing the local fields of the two spins. A strong local field causes the decoupling of the 
two spins, bringing them close to the product state, but the suppression of the purity is more dominant in two 
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FIG. 1: Schematic picture of the spin chains. We define the spins 1 and A'' as the focused spins and the spins which mediate 
the indirect interaction between the spins 1 and A'' as the media spins. We assume that we can tune the local fields on the 
spins 1 and N, while the other fields {/if }^"^^ are arbitrary but fixed. 

and three spin chains. On the other hand, in more than three spins, the decoupling by the local fields is more 
dominant than the suppression of the purity, and hence we cannot keep non-zero entanglement for any values 
of the local fields above a critical temperature. 

• The main factors of the entanglement enhancement due to the local fields are not only the increase of the purity 
but also their effect on the indirect interaction. To be more specific, the external fields affect the magnons 
which mediate the indirect interaction. The form of the maximizing local fields depends on this effect and has 
asymmetric forms in a particular parameter region. 

We show these results analytically and numerically. This paper consists of the following sections: in Section II, we 
state the problem specifically and give some definitions; in Section III, we give a general theorem on the entanglement 
enhancement which can be applied to any spin chains; in Section IV, we show the numerical and analytical results on 
the maximization of the entanglement in three-spin chains; in Section V, we show the numerical and analytical results 
on the maximization of the entanglement in four-spin chains; and in Section VI, discussion concludes the paper. 

II. STATEMENT OF THE PROBLEM 

First, we formulate the framework of the entanglement maximization problem and describe conditions. We consider 
a general XYZ A^-spin chain with external fields in the z direction. The most general form of the Hamiltonian of 
this system is given as follows: 

N-l N 

1=1 i=l 

where {(yl}i=x,y,z arc the Pauli matrices and we adopt the free boundary conditions. 

We hereafter consider the entanglement between the spins 1 and N at the ends of the chain. We refer to these two 
spins as the focused spins and to the external fields hf and on the two spins as the local fields. We refer to the 
other spins 2<i<iV— las the media spins (Fig. 1). The basic problem that we are going to solve is to maximize the 
thermal entanglement between the focused spins 1 and N by tuning the local fields /if and at a fixed temperature. 
We also fix all the other parameters, namely, { Jf , Jf , Jf } for 1 < z < — 1 and {/if } for 2 < i < iV — 1. We refer to 
the maximizing values of the local fields /if and as /iiop and /i^vop- 

Note that the maximizing local fields /iiop and /ijvop generally depend on the temperature T, or on the inverse 
temperature /3 = \/(kT) with k the Boltzmann constant. This is because we tune the local fields at a fixed temperature 
fj. Let us then define the high-temperature limit, in which we mostly develop the argument. In our high-temperature 
limit, we make /3||i?tot — h\al — /i^cr^H tend to zero, where || • • • || denotes the matrix norm. On the other hand, we 
let the maximizing local fields depend on /3 as we take the limit /3 — >■ 0. Hence, /3/iiop and PhNop may even diverge 
in our high-temperature limit. 

The density matrix of the total system in thermal equilibrium is 

Ptot = , (2) 

^tot 
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where Ztot = tr(e is the partition function. The density matrix of the focused spins 1 and A'' is 

Pin = triArptot, (3) 

where truv denotes the trace operation on the spins except the focused spins 1 and N. For the present system (1), 
the general form of the density matrix piN is given by: 



Pin = 



/Ptt F2\ 

pu Fi 

Fi 

\F2 pj 



(4) 



in the basis of the eigenstates of erf ® af^, where p^^, p^^, p^^, pn, Fi and F2 are real numbers. We have F2 = Q when 
Ji = Jf , ior 1 < i < N — 1, in particular. 

In order to quantify the entanglement, we here adopt the concurrence [18] as an entanglement measure. The 
concurrence C{pim) is defined as follows: 

C{pin) = max(Ai - A2 - A3 - A4, 0), (5) 
where {Xn}n=i ^'^^ the eigenvalues of the matrix 



^JpiN{cj\®al)plj,{cjl®cTl) (6) 

in the non-ascending order Ai > A2 > A3 > A4. For density matrices of the form (4), the concurrence C{pin) is 
reduced to the simpler form 

C{pin) = 2max(|i^i| - v^PttPII, 1^2 1 - VPtl^It'O)- (7) 
Then, the necessary and sufficient condition for the existence of the entanglement is given by 

max(|Fi| - VPttm> \F2\ - y/PnPit) > 0- (8) 

Thus, the present entanglement optimization problem for the spin pair (1,A^) is equivalent to finding the values of 
{/if , /i^} which maximize C{pin) with the parameters {Jf , Jf , Jf } {1 < i < N — 1), {h^} {2 < i < N — 1), and (3 all 
fixed. 



III. GENERAL THEOREMS ON ENTANGLEMENT GENERATION 



In the previous section, we formulated the entanglement maximization problem. In fact, there are cases in which we 
cannot generate the entanglement for any values of the local fields at all. In such cases, it is necessary to distinguish 
whether the entanglement is exactly equal to zero or rapidly approaches to zero. In the present section, we introduce 
a general theorem on a necessary condition to generate the entanglement by optimizing the local fields. In other 
words, we give a sufficient condition for the entanglement to be zero exactly for any values of the local fields. 

Theorem 1. Let us consider the XYZ chain (1) with {Jf , Jf , Jf } = { J^_i, J^„i, J^_i} = W , J^ , J''}- We tune 

the local fields hf and /i^, while the media fields {/if }^2^ ^^'^ arbitrary but fixed. There is a critical temperature 
above which the maximized entanglement between the focused spins 1 and A'^ is exactly zero if they are separated by 
two or more spins (N > 4). In other words, we cannot generate the entanglement for any values of the local fields 
above this temperature for TV > 4. 

Comments. This critical temperature yields a stronger restriction than the well-known ones [5] , which discuss the 
entanglement disappearance under a fixed Hamiltonian. Let us denote this critical temperature as Tc- We say Tc = 00 
if we can generate the entanglement at any temperatures, whereas we say Tc = if the maximized entanglement is 
equal to zero even in the ground state. 

Proof. In order to prove this theorem, it is enough to show 

max(Ff - PnPU' ~ PuPlt) < (9) 

after the maximization of the left-hand side with respect to hf and ft-^ in the high temperature limit /3 — > 0, where 
{-Fi, J2,Ptt'PT^iP^TiP4-i} elements of the density matrix defined in Eq. (4). Then, Eq. (7) yields the exactly 

zero concurrence in the limit /3 — 0. Since the system (1) has a finite number of degrees of freedom, the elements 
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{Fi, F2,p^^,p-\-l,Pl-\-,PH,} must be analytic as a function of (3. Therefore, there can be a finite value of /3 at which 
m.ayi(F^ — Pf-fP^^, F2 —p^Pif) = after the maximization. This gives the critical temperature. Note that the elements 
of Pin here are functions of /iiop(/3), hMop{P) and (3. 
Now, let us define the exponents ki and kat as 

/iiop = 0(/3-^^ ) and hNop = 0(;9-"" ) (10) 

in the limit /3 0, where ki and kat are real numbers. We estimate the order of each element of {-^1,-^2} and 
{P\tPii^P\iPi\\ following three cases: 

• Case (a): ki < 1 and kjv < 1- 

• Case (b): ki > kjv, ki > 1 and kjv > 0; or k^v > ki, ki > and kat > 1. 

• Case (c): ki > 1 and k^t < 0; or ki < and k^t > 1. 

Notice that the three cases cover the entire space of ki and kjv. 

Case (a). In this case, we can prove the inequality (9) by utilizing a necessary condition for the existence of the 
entanglement [19], which is 

trp?jv>^- (11) 

In the case (a), /3/iiop and Phfqop are of order f3^~'^^ and , respectively, and approach to zero in the high temper- 

ature limit /3 — >■ 0. The density matrix pi^ becomes proportional to the identity matrix because then /3||iJtot|| — > 0; 
hence we have 

limtrp2 =1. (12) 
^^0 4 

Therefore, in the casc(a), the entanglement between the spins 1 and N is exactly zero in the high-temperature limit. 

Case (b). To simplify the problem, we consider the case of /iiop,^Afop > 0, ki > fcjv, /^i > 1 and kn > 0, but we 
can prove the other cases in the same way. We define the exponent k as 

hiop - hNop = Oirn (13) 

in the limit /? — > and consider the case 

K > (14) 

in the following. We discuss the case k < in Appendix A. Here, in order to obtain the inequality (9), we prove the 
following; Fi and F2 arc of order of 

(9(^«i+«iv+«) and 0(;3''i+''"+'''), (15) 

or higher, respectively, where 

K ~ min(KAr, 1) > and k' — min(KAr, k, 1) > 0. (16) 

On the other hand, PttPj-J- ^^'^ PtlPlt ^^'^ both of order of 

C>(/32«i+2«"), (17) 

or lower. Then, the inequality (9) is satisfied below a certain value of /3. 

In order to prove (15) and (17), we separate the total Hamiltonian as follows: 

Htot = Hi + i^couplc + -f^modia + Hn, (18) 

where 

Hi —hio-pOl, Hn = ^AfopfATj 

TJ JX X X I TV I „Z I TX X X I jU„y „y I JZ Z Z 

Hcoup\e=J (Tia2 + J"(Jia:^ + J aia2 + J <Jn-i'^n + ^ '^n-i'^n + ■J ^n~i^n^ 

N-2 N-l 
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Now, we consider the term -ffcoupic, which couples the focused spins and the media spins, as perturbation. The 
unperturbed density matrix p[^l is given by 

pill = e-P"--P"«e-^"— (20) 

because -ffmodia and commute with each other. Because the external fields are applied in the z direction, the 
eigenstates of Hi + Hiq are given by {|titAf)i |tiiAf)i |iitAf)> IIiIa')} with the corresponding eigenvalues {— /iiop — 
/lATop, -hiop + liNop, hiop - /ijvop, hiop + hNop}; wc denote these eigenvalues as {eIJ^, Elj^, Efj^, We also give 

the eigenstates of i/mcdia as 

IV'Ldia) = SnlhMl^mN-l) +tnn2)Hl^)\iN-l) +UM^^^ (21) 

for n = 1, 2, • • • , 2^^^, where {IV'I^}: li^X^)' IV'i"'')} are the states of the spins from 3 to A^ — 2. Because the total 

Hamiltonian i/tot is a real matrix, the parameters {s„, t„, m„, w„} are real numbers. We define the corresponding 
eigenvalues of i/modia as {E^^^^^}. Then, the unperturbed eigenstates are given by 

{ItltTv) ® ICedia), ItllTv) ® IV'Ldia), \M n) ® ICcdia), IUIn) ® K.edia)}- (22) 

We define the perturbed eigenstates corresponding to each of (22) as {1^/;^?^): l''/'totf"'')j IV'tot'^)j l^tot^"'')}j respectively. 
We express them as 

ICf) = Itit^) ® ldia.n) + Itii^) ® ICidia.n) + \Mn) ® ICckit) + 1^1^^^) ® ICcku)' (23) 

for n = 1, 2, • • • , 2^-2 and ^ =tt, U, it, U, where | V'modia^tt) ' l^modia,ti) ' I V'mcdia^it) ^^'^ Idia.u) are the states of 
the spins from 2 to — 1 and may be not normalized. We also define the corresponding perturbed eigenvalues as 



{Kof , Kof , Kof , Kof }, which we express as 



Kof =£^f^+i?Ldia + '5Kof (24) 

for n = 1,2,-- - ,2^-2 a^j^^j ^ ^^^^ js^p^^ ^j^a^^ {-Bmcdial and {^Kofl are of order /3° because iJmcdia and 
-ffcoupie do not depend on the temperature /?. Then, we have 

^ = ^ = 0(/3°) (25) 

■^tot ^tot 

in the limit /? — > in the case ft-iop, /iwop > 0, where Ztot is the partition function of the total Hamiltonian. 

Now, we show the outline of the proof. First, we can calculate the elements {Fi, F2,Pff,p^^,Pif,pii} from (23) 
and (24) as 



^1 = ^ E E e"''''"- (Cfdia,tJ^n"lia4t)' (26) 



«=i {=tt,t4.,;t,U 
and 

m = ^E E ^-'''""'(CekttlCidia.n)' (27) 

*°* r, = l C=tt,ti,it,U 

for example. Second, we estimate the leading orders of {Fi,F2} and {ptJ-P^t'Ptt^'J-J-}- From the perturbation theory, 
we can obtain the approximate forms of {IV'tot)} and expand {Fi,F2} and {p'\-iPif,P'\-'fPii} with respect to /3. The 
elements {-^1,^2} are additive with respect to the indices n and ^; we define each contribution to the elements {i^i, F2} 
as 

^1 ~ \V^modia,til ^mcdia,it/' 

^2 = (V'media.ttl^media.u)' (2^) 
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Then, the elements {Fi, F2} are given by 

2iV-2 

^2 = ^E E e-^<'^V-«, (29) 
respectively. The elements {pfi-PitiPttPiii' ^'^^ other hand, are given by double summations as 



PtlPlt ^ ^2 ^ XV'mcdia.nl V^modia.tJ./ \ V'modia.it I %ncdia,it/ ' 



2 

i 

PnPu 



72 ^ XV-mcdia.ttl ^modia.tt/ \ V^mcdia.iJ. I V^modia,4.4,/ ' 



Because the summands of (30) arc all positive, we can obtain the following inequalities for {p^^\.p\.■\^p^^^^p\^\^'■ 

PTWT - 72 ^ \V^media,t4,l 'f^media.tJ./'^ \ rmedia,^t I ^media,^t/ ~ t4-,4.t' 



tot 



n=l 



PnPii - y2 ^ \V^media,ttl ^media.tt/ \ ^media.U I ^modia.U/ ~ tt,U ' ^ ' 

^tot 

where we pick up only the terms oi n ~ n' and = =tt from the summations in (30). The elements ^-^^ and 
■^tl a are additive with respect to 72; wc define each contribution to the elements {P^^^^i ^-f/^^} 

t-l,4.t ~ \''^modia,t^| V^modia,t^/\'f^modia,itl ^mcdia,^t/' 

tt4-l ~ \''^mcdia,tTl ^media,tT/ \'f^media,4,4,| V^media.U/' ^ ' 

Then, the elements {P^J ^^^^P"^^ ^;^} are given by 



N-2 



1 2 

■f^' '* - 72 Z^ "2 -f^ 



2Af-2 

tt,U ~ 72 2^ ^ -^tt,U' 



^to 

respectively. In this way, we calculate each contribution to the elements {^1,^2} and {Pj'^j^^, P^^^^j^} separately. 

In accordance with the above outline, wc first calculate jV'tot'^^), namely the perturbed state of Ititw) ® IV'media)- 
In order to simplify the calculation, we rewrite the perturbation Hamiltonian -ffcoupic as follows: 

iJeoupic = + l)olal + (1 - l)o\ol\ + ,ra\al 



where 



: J{a\a^ +(7^0-++ -i{at^t + ^\ ^2 )] + J'' ^1^2 

+ J{a%_^a^ + cr^_iO-+ + 7((T+_^(T+ + o^_^a^)\ + JV^_iCr^, (34) 



7^1: IV 

j^r + jy, 7 = - — -. (35) 

' jx JV ^ I 
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From the calculation in Ref. [20], the leading terms of the elements {IV'mcdil?)} ^ =11; tl. It; 11 in Eq. (23) are 
given by 

IC%ltt) = ICcdia), (36) 

IC"edVa,Ti) =^^(75nolt2)|^J;I)|;A.-l) +t„„|t2)|^,ti)|t^-l) 
—ZhNop^ 

+ 7"nol;2)|^iI)|;iV-l) +^«„ol;2)|^it)|tiV-l) +0(r~)), (37) 

iCedlit) =r^(7^"oi;2)k'i;i)itiv_i) +7in„i;2)iv;j;i)i;iv-i) 

+ «no It2) k'il) ItiV-l) + \iN-l) + OiP'^- )) (38) 

and 

i&lu) =77-^(7'^noii2)i^:i)i;^-i) +7i«oi;2)i^:i)it^-i) 

4/7,iopn.Arop \ 

+ 7"n„lt2)|^il)|;Ar-l) +^i'„„|t2)|V;it)ltiV-l) +0(r~)). (39) 

Now, we calculate the contributions of IV'tot^^) to the elements {i^i,F2} and {P^^^ ^^,P^^ which are defined as 
|^«o,tt^^no.tt| {Pl^^^J^^P^^'^ll}. From Eqs. (36)-(39), we have the elements {FI'°^^\f^"-^^} as 

iriiopllNop ^ 

+ 7'<no"no(^iIl'^It> + IWr.o^noii'llli'ii) + 0(/3^")) (40) 

and 



"o,tt 



^ 4/liop/lArop 

+ 7i«o«no(^:il^i!) +0(r")) (41) 



as weU as the elements {P^^"J_^ , P^^'J^} in the forms 

P. 



«+/?-<5P;-'^) (42) 



-'■"J'nop''-Arop 



and 



t4 

p"o,tt 



^"'nop"Arop 



(7'4o + iXo + + < + )) 



- ^< + '^^"t.U ) (43) 

^""lop"-Arop 

in the limit /3 — )• 0, where we define the real numbers of order /3° as 5P^'^'^^ and so as to satisfy the above 

inequalities. 

Then, we sum the elements {F'^'^^ ,f:^'^^} and {P"J/^, P"^^^} with the Boh zmann weight e f'^tot over the label 
n. First, we calculate the summation of each of {P"'^^, P"'^^}. Because the spins 2 and [N — 1) are separated by 
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{N — 4) spins, the correlation between the spins 2 and (A'' — 1) are generated by the {N — 3)th-order perturbation of 
i7„iodia- Therefore, we obtain 

K'^^-i>o = 0(/3"^), (44) 

where (■ ■ ■ )o denotes the thermal average with respect to p['^^ in (20) and ai > A'' — 3, a2 > A'' — 3. Since we are 
considering the case iV > 4, we have ai > 1 and a2 > 1- This is the key to the fact that the entanglement can vanish 
for > 4. From the equations 



(o-^tT^_l + cr^_i)o 



2 

- r| e 2 



=tr(e-'^^'— |;2tiV-i)(t2;iV-i|) 



— E e-^^"--(t2;^-l |V'S.cdia)(«edia|;2tAr-l) (45) 



and 



4 
2 

_, I -/3g„wia li2lAr-l)(t2tAf-l| + |t2tjV-l) (l2lAf-l | 

-trie 2 



=tr(e-'3^— |t2tiv-i)(;2;w-i|) 

2^-2 

7^ E e-^^»--(;2;Ar-l |<.cdia)(Ccdia|t2tiV-l), (46) 
'^mcdia 



n=l 



where 2'mcdia = tr(e /^-ffmcdia^^ (j+ = (^^x _|_ i(jy)/2 and <t = (cr^ — icr^)/2, we also have 

E e-^^"--t„^.„(^1;^|^f ) - Oin (47) 



and 

2 



N -2 



E e-^^»--s„u;„(^i^|V;r> = 0(r), (48) 



n=l 



where a = min(ai,a2). Moreover, because {5El^^^} in Eq. (24) are of order we have 

e-/5sro'^ ^ e-^(^'«+^."-»' (l - /3<5i?t"f + 0{P^)) , (49) 
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forn = 1,2, - •• ,2^^2 ^ rj.-^^^^ obtain from Eqs. (40), (41), and (47)-(49), 

2iV-2 

2iV-2 2-^-2 

n— 1 n— 1 

^e-^^'^- (0(/3") + 0(/?^« )) + e--^^^" -— 1- 0(/3) 

4/liop/lArop ^ ^ 4ftiopn.Arop 

=e-/'^^S X 0(^'^i+'^"+«) (50) 



and 



2JV — ^ 

Y e-^^"-^^^' F^^^^ = e-'3^^«— — f0(/3") + ©(/J'^") + 0(/3)) = e"''-^!" x 0(^'=i+«iv+«)^ (51) 

4/llnnAlAfnn V / 

n— 1 



4/liop/lJVop 

where a = niin(ai,a2) > 1 and k is defined in (16). 

We can similarly calculate the contributions of the other states {|^/'"of)} =tiiitjii)- From the calculation in 
Ref. [20], the leading terms of the elements {IV'me'dia 4)^ ^ =ttjtiiit7ii in Eq. (23) can be obtained, and we have 
the contributions of the states {iV'tot^"'')} to Fi and F2 as 

2^-2 

y e-^^"- V;'^^ = e-^^^« — fo(/3") + 0(/3«") + 0(r) + 0(I3)) = e-''^i'« x 0(/?«i+""+"'), 



4/liop/lArop 

y e-^K^^'-F^'ti ^ ^-PEl^J^ (0{ir) + 0(/3''") + 0(/3)) = e'^^^" x 0(/3"i+«"+«), (52) 



^^-^ 4/liop/lArop 

where we utilized Eqs. (47)-(49) and k' is defined in (16). From the inequality 



Ztot 

for ^ =tt,ti,it,ii, we finally obtain (15) by substituting Eqs. (50)~(52) into Eq. (29). 

□n^tt pn-tt 



< 1 (53) 



Second, we calculate the summation of each of {i^IV ^t-, -Pft ii}- "^^^ parameter Wn in Eq. (21) cannot vanish for all 



n. Therefore, we have 



n=l 

-,iV-2 



^tot 

^E^"''''""^" = O(/30), (54) 
where we utilized Eq. (25). Hence we obtain 

2N"-2 

> / --'^^""K1 + /5''"^n"i4t) = 0(/3-^+-«) (55) 



and 



> ,J.. E e-^^^""K + r-'^^n-'u) = 0(/3-^+-«). (56) 

-^"'nOp'^ATop ^tot 
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Because p^^Pi^ > PjJ a-nd p^i^Pn > P^^ wc obtain (17). Thus, the inequahty (9) is satisfied below a certain 
value of /3 in the case (b). 

Case (c). To simplify the problem, we consider the case of /ii > 0, ki > 1 and kat < 0, but we can prove the other 
cases in the same way. In this case, we prove the following; and are both of order of 

0{I3^'''+^), (57) 

or higher. On the other hand, PffPn and p-fiPif are both of order of 

0(/32-i), (58) 

or lower. Then, the inequality (9) is satisfied below a certain value of /3. 

In a similar manner to the case (b), we separate the total Hamiltonian as follows: 

^^tot = Hi + -ffcouplc + -f^modiai (59) 

where 

N-l N-1 
^Ccdia = E (J^'^f'^f+l + J1«1 + Jt^t^Ul) + E ^< + (60) 



Note that in the case (c), the norm of the Hamiltonian H^^^^^^ is of order /3° because kjv < 0. Here, we consider 
the interaction term ^couple which couples the spins 1 and the other spins, as perturbation. Then, the unperturbed 

density matrix p^^'J' is expressed as 

^e-^^^e-'^^"-, (61) 



The eigenstates of Hi can be given by {|ti), with the corresponding eigenvalues {— /iiop, ^lop}; we denote these 
eigenvalues as {eI,Ei}. We also define the eigenstates of H'^^^^^-^ as 

ICedia) = <lt2)|0r>lt^> +4lt2>|^I^>|;Ar> +<|;2>|0f )|t^) +<|;2)|^t^)|;A.) (62) 



for n — 1, 2, • • • 2^ ^, where {l^ij^^), \<j>li^), l^n"'')} are the states of the spins from 3 to — 1. Because the total 

Hamiltonian Htot is a real matrix, the parameters {s[^,t'^,u'^,w'^} are real numbers. We define the corresponding 
eigenvalues of ^^^.^dia as {£^,nodia}- Then, the unperturbed eigenstates are given by {|ti) ® l^modia)} and {||i) (g) 
I'^media)}- define the corresponding perturbed eigenstates as {|0"otf)} and {|(/)"ot^)} and the corresponding perturbed 
eigenvalues as {E"^^ , El^^^} . Wc express them in the forms 

IC?) = Itl) ^ l&a,t) + lil) ® ICedia,^) (63) 

and 

<r = i?i'+i?:"edia+'5i?;r (64) 

for 71 = 1, 2, • • • 2^^^ and 77 =t, I, where l^^'o^ia f) and \(l)l'^cdia i) are the states of the spins from 2 to N and may not 

be normalized. Note that {E^^^^-^} and {SE^^f} are of order /3° because ||-ffmcdiall and H-ffcoupidl are of order /3°. 
Then, in the limit /? — > 0, we have 



e 



= 0(/3°) (65) 

^tot ^tot 

in the case /iiop > 0. 

Next, we calculate the elements {^1,^2} and {pti-P^t'PttPul- The calculation in [21] gives the state |0"ot'^) the 
form 

\&) = Itl) ® ICedia) + ^5^iii) ® (74„i;2)i0];i>itAr) + 7Ci;2)i<^:i;>i;A.) 

+ <olt2)|^iI)lt^) + <lt2)|0it)i;^) + o(r^)). (66) 
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Then, we can calculate the contribution of |^"ot'^) to the elements {Fi,F2} and {ptt'-PtJ-'-P^T'PJ-l}; denote them 
as {Ff°'^F2"°'^} and {p'^f^p'^f^pl^'^^plf}. First, the elements {Fl""'^ , F^"^^} arc given by 



1 -z^(74o<('/'itl'^n!) + + (67) 

and 

{<o<o i<t>il + {<t>ll\^ii) + )) . (68) 

"o,t "o.t "o.t "o.tl 



Second, the elements {p""^ '^'ti '^"t '^'u J ^^"^ given by 



tt ^ '^0 ' no ' 
no,t j-'2 I '2 



^"■lop 



72 

no,t 



-i^(7'C+<+0(/3-)). (69) 

Then, we sum and with the Boltzmann weight 6^'^^'°*^ over the label n. First, 

we calculate the summation of each of {-F"'^, -^2''^}- From the same discussion as in Eqs. (47) and (48) in the case (b). 



we have 



and 



2 



N~l 



^ e-^^»-»C<(0ti|^it) = 0(/3"') (70) 



n=l 



2 



N~l 



E e-^^.-.^s>;(^i^|0r) = O(r'), (71) 

ri=l 

where the exponent a' is defined as follows: 

a = mui(ai, aj), 

K^^) = O(r'0, K^^) - 0(/5"^), (72) 

where a'l > N — 2 and a'j > iV — 2. Since we are considering the case > 4, we have a'l > 2 and a'j > 2. On the 
other hand, because {SE^^^^} in Eq. (64) are of order /3°, we have 

e-^^^^L"'" ^ e-'3(^;'+<"-.) (1 - /3,5ii;Jt'' + 0{f3')) , (73) 

forn = 1, 2, • • • , 2^-1 and =t,i. Then, we obtain from Eqs. (67), (68), and (70)-(73), 

2iV-l 



n=l 



n— 1 n— 1 



— 2/iiop V / — 2ft,iop 



and 



n— 1 ^ 
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We similarly calculate the contributions of the other states {iV'tot^)}; then, we finally arrive at (57), 

F2 = ^ e-'^^i X 0(ri+i) 0(/3'=^+i), (76) 



■^tot . , 

where we utilized the inequality 



< 1 (77) 

■^tot 

for ?7 =t,;. 

Second, we calculate the summation of each of From Eq. (69), we obtain 

2«-l 

Y: e-^^-'V;f = e-^^^ X 0(/3°), 

Tt = l 

2iV-l 

n=l 

2iV-l 

J2 e-^'^'-Vlt^ = ^^^""^^ ^ 0(/32«i), 

n=l 

2iV-l 

^ e-^^-^p^^ ^ X 0(/32«i ), (78) 



n=l 

because sj^ + m^, + w^^, 7^s^ + and 7^t^ + cannot vanish for all n. From Eq. (65), we have 



'^tot — ; 

n— 1 

2iv-i 

m>^E^"''^^""prt' = o(/5"^)' 

^tot — ; 

n— 1 

.-,7V-1 



> ^ E ^-''''^""nY = )■ (79) 

We thereby obtain (57). We thus obtain (57) and (58) and hence the inequality (9) is satisfied below a certain value 
of /3 in the case (c). 

Thus, we prove the inequality (9) in the cases (a), (b) and (c). This completes the proof of Theorem 1. 

In this theorem, we discuss the case of {Jf, Jf, Jf} = {J^-n ^n-i^ Jn-i) = {J^ , , J^}, but the case of 
{Jf, Jf, Jf} 7^ {J^_i, J%_i, J^^i} can be also proved by extending this following proof; we only have to repeat 
the same calculation by changing -ffcoupic in Eq. (19) into 

i^couplc = Jl(Tf(T2 + Jl^l^^l + Jl^^l^l + Jn-1'^N~1'^% + Jn-1<^1-1<^1 + Jn -I'^N -I'^N ■ (80) 

IV. ENTANGLEMENT MAXIMIZATION IN THREE-SPIN CHAINS 



In the present section, we discuss the maximization problem in three-spin chains. A significant result of the present 
section is that the maximizing local fields are asymmetric to each other in a particular parameter region. This effect is 
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FIG. 2: Schematic picture of the three-spin chain. We define the spins 1 and 3 as the focused spins and the spin 2 as the media 
spin. The spins 1 and 3 are symmetric to each other and they indirectly interact with each other through the spin 2. 



nontrivial; the positions of the spins 1 and 3 are geometrically symmetric to each other if { Jf , Jj*, } = {J2,J2,J2} 
and the entanglement is invariant with respect to the exchange of the two spins. Nevertheless, we obtain |/iiop| |/i3op| 
in a region. We show that the asymmetry can be mainly attributed to the behavior of the magnon. This effect indicates 
that the purity of the focused spins is not the only criterion for the enhancement of the bipartite entanglement. This 
is different from the case of the entanglement maximization in two-spin systems [14]. 



A. Numerical results 



The Hamiltonian which we consider in the present section is 

2 
i=l 

+ hlal + K,,^,^CT^2 + hl^l (81) 

where the spins 1 and 3 are the focused spins and the spin 2 is the media spin (Fig. 2). We assume > > J^. 
We solve the entanglement maximization problem about the spins 1 and 3 by fixing the temperature T and the field 
^modia ^^'^ Spin 2. In order to solve this maximization problem numerically, we use the random search method and 
the Newton method together. According to Theorem 1, if the two spins were separated by two spins, there would 
always be a critical temperature above which the maximized entanglement is exactly equal to zero. In the present 
case, the focused two spins are separated by only one spin and therefore Theorem 1 docs not apply; the elements 
1^1 1 7 y/PttPil-' 1^2! and ^/ptJpJt are shown to be of the same order in the same way as in the proof of Theorem 1. 
Therefore, it generally depends on the interaction Hamiltonian and the positions of the focused spins whether the 
critical temperature exists or not. As for the Hamiltonian (81), except the case of the Ising chain with /i^jg^ia = 0, we 
can prove [22] that the entanglement between the spins 1 and 3 can exist at any temperatures by letting /if = /i| — >■ 00. 

In Fig. 3, we show the numerical results on the entanglement maximization. The main feature is that in a certain 
parameter region the maximizing local fields /iiop and /isop are asymmetric to each other, namely jft-iopj 7^ |^3op|- 
In this region, the asymmetry must be essential to the enhancement of the entanglement. The asymmetry appears 
continuously (solid line in Fig. 3) or discontinuously (broken line in Fig. 3). Note that in these cases there is no 
critical temperature above which the maximized entanglement would be zero. 



B. Analytical argument 

Here, we argue the origin of the asymmetry for the XX model; the phase diagram of the other models are not 
much different from the XX model qualitatively. The main reason of the asymmetry is the strong dependence of the 
indirect interaction on the local fields. Then, we focus on the magnons which mediate the indirect interaction and 
discuss the effects of the local fields on the magnons. We show that the following three points affect the asymmetry: 

1. The Boltzmann weights of the states with the magnons. 

2. Localization of the magnons. 

3. Suppression of the off-diagonal elements {-Fi, F2} due to mixing of different magnon states. 
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(c) 



FIG. 3: The parameter ranges in which the asymmetry appears (a) in the case of J"^ = = 1 and (b) in the case of 
— 1,J^ — 0.5. In the 'Asymmetric' phase, the optimizing fields /iiop and /isop are asymmetric to each other as |/iiop| 7^ I^SopI, 
while they satisfy |/iiop| = |/i3op| in the 'Symmetric' phase. On the solid line, the asymmetry appears continuously, while on 
the broken line the asymmetry appears discontinuously. The ratio /iiop//i3op is shown for (c) ftmedia = 4 and (d) T = 1.5 on 
the chained lines of (a). 



We define a magnon as a spin flip; a down spin in the background of up spins or an up spin in the background of 
down spins. For example, the magnon number is two for both of the states Itttii) and jftiii)- In the XX three-spin 
chain, the number of the magnons can be either zero or one according to the above definition. The magnon does not 
exist in the states |tit2t3) and |ii4.2i3), nor the entanglement. This suggests that the magnons mediate the indirect 
interaction and is essential to the entanglement between the focused spins. Even if the magnons exist, however, the 
entanglement can be very small when the magnon is localized in one site. It is also possible that the off-diagonal 
elements {Fi,F2} in the density matrix (4) are suppressed by the mixture of the two kinds of magnon states. As 
shown in Eqs. (45) and (46), these elements are related to the correlation between the focused spins; 



^1 = 



4 



^2 = ^ • [aZ] 

Let us show an example; each of the two states 

^(Itit2;3) + i;it2t3)), ^(itit2;3) - i;it2t3)), m 

has one magnon and its entanglement is maximum. However, the correlations which their magnons mediate are 
opposite to each other; the first state has (erf erf ) = (crfcrf) = 1, whereas the second one has (crfcrf) = (cicrf) = — 1. 
Then, the off-diagonal elements {-^1,^2} are also opposite to each other. Therefore, the entanglement is enhanced if 
we suppress the mixedness of the one state of (83). 
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Let us analytically show how the above three points affect the asymmetry. First, we consider the case of T ~ and 
^mcdia ^ that is, the upper left corner of the phase diagram in Fig. 3. In this case, we have /iiop = h^op > in the 
'Symmetric' phase and /iiop, ^3op > but /iiop 7^ /i3op in the 'Asymmetric' phase. Let us consider the interaction J 
and the local fields /iiop, h^op on the focused spins as perturbation. Then, the four unperturbed ground states |tit2t3)j 
|tit24-3), Iiit2t3) and |iit2i3) are degenerate; we define these four states as Because we assume /i^edia ^ ^ ^ 

we consider the mixing of only these four states. The first-order contribution of the excited states |tii2t3)) |ti4-2i3), 
Iiii2t3) and 14,14,2^3) to the above four unperturbed ground states are given in the forms 

I'/'l) = Itlt2t3), ei = -hlop - h^op - /iLdia: 

|^2> = /C+(|tlt2;3> + j|ti;2t3)) +/C_(|;it2t3> + j|ti;2t3)), 



£2 = -(5 - \ {hlop - hsop)^ + - 



media' 

1^-3) - /C_(|tlt2;3> + jltli2t3)) - IC+ilUhh) + jlti;2t3)), 
£3 = + \/ (hlop - /l3op)2 + S'^ - ft.n,edia> 

ji;i;2t3> + ji 



li^i) - i;it2;3> + -i;i;2t3> + -iti;2;3>, £4 = hi^p + ha^p - Ccdia + 2(S, (84) 



where 6 = - JV(2^modia) and 



IC^ = ±. l±—Jl^^J^^^=. (85) 
\/2Y Vihiop-h,opr+S^ 

We thereby calculate the elements Fi, p-^-^ and pn of the density matrix pi2, which contribute to the concurrence. 
First, the leading terms of and pn are obtained as follows: 

Pn = — and = -7} • (86) 

^tot -^tot 

Therefore, ^Jp\\Pii is equal to e^^^"^'"^'^'^~^^ jZ tot and does not depend on h\op and hzop- Then the entanglement 
depends on h\op and h^op mainly through F\ as is shown in (7). 
We obtain the leading term of F\ as follows: 

Fx = /C+/C_ 

^tot 

Let us see how Fi depends on the asymmetry /iiop — h^op- We consider the following three terms of Fi separately: 

{h,op-h,opr 



2 (hlop - /i3op)2 + 52 



We start from the point /iiop = h^op and gradually increase the asymmetry hiop ~ ^3op to see how the above three 
terms affect the entanglement through the element Fi. We argue that there are two effects competing with each other, 
which may yield the phase boundary in the upper left corner of the phase diagram in Fig. 3. 

The two factors /C_|_ and /C_ on the right-hand side of the first term in (88) appear as coefficients in the states 1-02) 
and \ip3) in (84). As the asymmetry hiop — ft.3op is increased in Eq. (84), the factor /C_ decreases and therefore the 
states \tp2) and 1-03) become dominated by the states |tit2i3) and |iit2t3): respectively; the factor /C+ decreases when 
the asymmetry is reversed to h^op — hiop. Then the magnon becomes localized on the site 1 or 3 as the asymmetry 
hlop — hsop is increased. This is a negative effect on the entanglement; the entanglement would be suppressed if the 
magnon does not move around the system. Indeed, the element Fi in (87) and hence the entanglement in (7) can 
decrease because of the decrease of the first term in (88). 

The second and third terms in (88) arc the Boltzmann weights of the states \ip2) and |'03)- As the asymmetry 
hlop — h^op is increased, the second term increases while the third term decreases exponentially. The state \ip2) has 
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positive transverse correlations {ufa^) and {a\a^), while the state |-03) has negative ones. As has been shown in the 
example of (83), the increase and decrease of the respective Boltzmann weights enhance the transverse correlations. 
This is a positive effect on the entanglement; the transverse correlations can promote the off-diagonal elements {^1,^2} 
in Eq. (82). Indeed, the element Fi in (87) can increase because of the exponential increase and the decrease of the 
second and third terms in (88). 

The above negative and positive effects of the asymmetry on Fi and the entanglement compete with each other. In 
other words, in some cases the entanglement is maximized without the asymmetry, while in other cases it is maximized 
by introducing the asymmetry. This may be the reason for the phase boundary in the upper left corner of the phase 
diagram in Fig. 3. 

Second, we consider the case /imcdia > T 3> 1, that is, the phase of the asymmetry in the upper right area of the 
phase diagram in Fig. 3. We numerically obtained hi^p ^ J, h^op 'ti — J and |/iiop| 7^ |^3op| in this area. We show 
that the asymmetry indeed promotes the entanglement by starting from the point ft-iop = — /isop and increasing the 
asymmetry hi^p + ft-sop gradually. We argue that the increase of the entanglement is due to delocalization of the 
magnon. In order to calculate the entanglement in the present case, we regard the interaction between the focused 
spins and the media spins as perturbation. 

Because |ft-iop| ^ 1 and |ft.3op| ^ 1, we consider only the two states |tit2i3) and |tii2i3) as the unperturbed states. 
Then, these two states in the first-order perturbation are given by 

J 

\lpl) = Itlt2l3) H TUTT 1 I ,z 7ltll2t3> + T7TT 1 . ,z TTT " FT rT|llt2t3), 

-2(|/l3op| + /l^^diJ 4(|/l3op| + /iLdia)(^lop + l^3op|) 

= ~'^lop - |^3op| - ^mcdia (89) 

and 

J J2 

|V^2) = ItlMs) + -^777 — TI 7lllt2l3) + T7T — wT — FT |Tlill2t3), 

-2(/llop + Vcdiai 4(/liop -I- Vcdiail/^loP + I "'Sop | ) 

£2 = -^lop - |/l3op| + /imcdia- (90) 

From these expressions, we can calculate ^p^tPii 

^;p^^,^(/..op + |/.3o.|) J' (91) 

4U'l3op| + "-media) I "lop + Vodiai 

Similarly, the element Fi is given by 

\^4(|/i3op|+/lLdia)(/ilop + |/i3op|) 4(/liop + /ir„edia)(/»lop + |/»3op|) )' ^ ^ 

To clarify the effect of the asymmetry /iiop + h^^p, we put 

^lop ^ ho + 5h and h^op = -/^o + Sh, (93) 
where ho ^ 5h > 0. Then, Eqs. (91) and (92) arc recast into 

4L(/lmcdia + ^0)^-<^/l^J 

and 

p ^ g2/3/io [ ^! e'3'''-d» + — e-'3''.-d.a 1 f95) 

\Sho{ho + h;^^^,,^~5h) 8ho{ho + h:^^,,^ + Sh) J' ^ ' 

The former ^/pttPii does not depend on the symmetry Sh up to the first order. On the other hand, the latter Fi and 
hence the entanglement in (7) increase as Sh is increased; the increase of the first term on the right-hand side of (95) 
excels the decrease of the second term because of the difference in the Boltzmann weights. 

The increase of Fi can be related to delocalization of the magnon. Since |/iiop| 3> 1 and |/iArop| S> 1, the magnon 
is almost localized in the site 3 in \tpi) of (89) and in the site 1 in |'(/'2) of (90), but slightly delocalized owing to the 
perturbation terms. As we increase Sh, the perturbation terms of increase, but those of \tjj2) decrease. Therefore, 
the magnon is delocalized more in IV'i), while it is localized more in \ip2)- This corresponds to the increase of the first 
term and the decrease of the second term in (95). Because the Boltzmann weight of the state j^'i) is greater than that 
of |'02), the magnon delocalization in \ipi) excels the magnon localization in \ip2)- This corresponds to the increase of 
Fi. 
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FIG. 4: Schematic picture of four-spin chains. We define the spins 1 and 4 as the focused spins whereas the spins 2 and 3 
as the media spins. The spins 1 and 4 indirectly interact with each other through the spins 2 and 3. In this system, the two 
focused spins are separated by two spins. Therefore, according to Theorem 1, there is a critical temperature above which the 
maximized entanglement is equal to zero. 
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FIG. 5: The phase diagrams for the four-spin chain, (a) in the case of = = 1 and (b) in the case of — 1, = 0.5. 
In the 'Asymmetric' phase, the maximizing fields /iiop and /i4op are asymmetric to each other as |feiop| 7^ |/i4op|, while they 
satisfy |feiop| ~ |^4op| in the 'Symmetric' phase. The asymmetry appears continuously on the solid line, while it appears 
discontinuously on the broken line. The maximized entanglement vanishes beyond the chained line. 



V. ENTANGLEMENT MAXIMIZATION IN FOUR-SPIN CHAINS 



In the present section, we consider the maximization problem in four-spin chains (Fig. 4). As has been proved 
in Theorem 1, the four-spin chain is the shortest one in which the end-to-end entanglement cannot be generated in 
the high-temperature limit. We mainly discuss the critical temperature and its dependence on the interaction of the 
spins. 



A. Numerical results 



In the present section, wc consider the XY spin chains given by the Hamiltonian 

3 

J/tot=5^(JVfa,?+i + JVfa,:Vi) 
1=1 

+ hlal + hlal + /i^,edia('^l + ^D- (96) 

We solve the entanglement maximization problem about the focused spins 1 and 4 by fixing the temperature T and 
the field /i^edia '^'^ media spins 2 and 3. In order to solve this maximization problem numerically, we used the 
random search method and the Newton method together. According to Theorem 1, there always exists a critical 
temperature above which the maximized entanglement is exactly equal to zero because the focused spins 1 and 4 are 
separated by two spins. 
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We show the phase diagram of the XX spin chain and the XY spin chain in Fig. 5. As in the case of the three 
spins, there is an 'Asymmetric' phase where the maximizing local field |/iiop| is not equal to |/i4op|- Moreover, the 
maximized entanglement is equal to zero in a region. The qualitatively different behavior of the critical temperature 
between the XX and the XY chains is due to the conservation of the angular momentum in the z direction, as we 
will argue in Section VB. 

B. Difference between the XX and the XY model 

We discuss the behavior of the critical temperature in the XX and XY chains. In the XX chain, the critical 
temperature increases as the media field /i^edia increased, while in the XY chain it does not. This difference is 
attributed to the conservation of the angular momentum in the z direction. In the XX chain, we can suppress the 
mixture of the states with more than two magnons by choosing the fields as 

hi = -/lo 

hi = ho, for i = 2,3, •• - iV, (97) 

and ft-o/? ^ 1- Then, the density matrix e~'^^'°' is almost equivalent to the ground state of Htot and the mixture of 
the other states can be suppressed exponentially by increasing h^. The ground state is given by the following form; 

aollit2 ■ ■ • tjv) + ai|til2 • ■ ■ tw) + a2|tit2l3 ■ ■ • tjv) 

+ ---+aiv-i|tit2 •••Ijv), (98) 
where we can calculate by the fcth order perturbation to have 

a.^og) (99) 

for fc = 1, 2, • • • iV — 1 with the factor h^ coming from the energy denominator. In the ground state (98), the element 
Pa is equal to zero because there is no state with more than one down spins in (98). The entanglement between the 
spins 1 and N exists because cx |ajv-i| oc |J//io|^~"^ > 0. The mixture of the excited states generally destroys 
the entanglement, but is suppressed exponentially because of the Boltzmann weights. The entanglement between the 
spins 1 and N thereby survives. 

On the other hand, in the XY chains, we cannot control the number of the magnons in the ground state by 
increasing the media fields. Therefore, PttPii is not zero in the ground state, which invalidates the argument for the 
XX model. This may account for the fact that the critical temperature does not increase as the media fields /i^jcdia 
are increased. 



C. Calculation of the critical temperature 

Let us calculate the critical temperature analytically. It is generally difficult to solve the entanglement maximization 
analytically. At the critical temperature, however, we numerically confirmed that the maximizing local fields take 
the form /iiop = /i4op — ^ oo in some region of the 'Symmetric' phases in Fig. 5. For the XX model, this happens 
everywhere in the 'Symmetric' phase in Fig. 5 (a), whereas for the XY model it happens only in the lower one of 
the two 'Symmetric' phases Fig. 5 (b). In such regions, we can derive the critical temperature analytically. In the 
'Asymmetric' phases in Fig. 5, on the other hand, the maximizing local fields do not have simple forms. 

After the calculation in Ref. [23], we obtain a necessary and sufficient condition (8) for the existence of the entan- 
glement in the form of the following two inequalities: 

- mm) = [45_5+7sinh/?|e3| + (1 + 7') sinh/3|ei|]' 
-Ztot[(7'52 )e-^^^^ + (52 +7*5^)6-^^^ +272 cosh/?|ei|] > 0, (100) 

or 



ZtoAFi ~ PuPi^) = [25_5+(l +72)sinh/3|e3| +27sinh/3|ei|] 

- [{I'^Sl + SDe-f"'' + (5! + 7^52 )e-^^^ + (7^ + 1) cosh/?|ei|] > 0. 



(101) 
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FIG. 6: Comparison of the phase boundary (a) in the case of = J'' = 1 and (b) in the case of = 1, — 0.5. The critical 
temperature calculated according to the inequality (100) is indicated by the chained lines and the numerically calculated one 
by the solid lines, which are the same as in Fig. 5. The critical temperature calculated by the inequality (101) is also indicated 
by the broken line in (b). In the XX spin chain, the approximation gives almost the same values as the numerical ones in the 
'symmetric' phase, while in the XY spin chain, the approximated values from the inequality (100) fit the numerical ones in the 
range of < ftmcdia < 0.851. 



where J and 7 are defined in (35), 




2h^ 

'""'^'^ (102) 
and 

ei = -62 = -J" - J^ 

Ztot = e-^'' + e-P'^ + e-P'^ + e-P'\ (103) 

In the case of the XX model with 7 = 0, in particular, the above condition becomes simpler; because we have 
F2 = for 7 = 0, the inequality (101) is never satisfied. The condition (100), on the other hand, reduces to 

sinh2(/3|el|) - Ztote-^'^ = sinh2(2/3J) - Ztote""'^''-^- > 0, (104) 
■1- sinh^(2/3J) 



or 



2cosh(2/3J) +2cosh(2/3/i 



> 1. (105) 



media/ 



In Fig. 6, we compare the critical temperatures calculated according to the inequalities. (100) and (101) with the 
numerically calculated ones. The approximations generally give precise estimates for the XX model in the 'Symmetric' 
phase. For the XY model, the approximate estimate shown in Fig. 6 (b) is due to (100) in the range bounded by the 
broken line and is due to (101) in the range bounded by the thin solid line. The former gives good estimates in the 
range of < /i^.^;, < 0.851. 



VI. SUMMARY AND CONCLUSION 



We have analytically and numerically studied the maximum value of the thermal entanglement between two spins 
which indirectly interact with each other. We showed the general theorem on its temperature dependence. We proved 
that the maximized entanglement is equal to zero above the critical temperature if the two spins are separated by two 
or more spins. This result tells us that as the temperature goes higher the thermal entanglement between distantly 
separated two spins cannot survive with any local controls. 
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Secondly, in the three-spin chains and the four-spin chains, we showed properties of the maximized entanglement 
which is calculated numerically. In the three-spin chains, we showed that the maximizing local fields are not symmetric 
as /iiop 7^ h^op (Fig. 3) in some parameter regions. We have explained the asymmetry qualitatively and quantitatively. 
We attributed the appearance of the asymmetry to the effects on the magnons, which mediate the indirect interaction. 
In other words, the local fields affect the effective interaction between the focused spins, while the direct interaction 
is determined independently of the external fields. In the four-spin chains, we demonstrated that the maximized 
entanglement vanishes above the certain critical temperature. Because of the difference in the symmetric properties 
between the XX and XY spin chains, the dependence of the critical temperatures on the media fields are qualitatively 
different between the two systems. We calculated the critical temperature analytically in some parameter regions. 

In conclusion, we have clarified several properties of the entanglement which is generated from the indirect inter- 
action. Our study has given one of the general limits for the entanglement generation. We have also shown some 
properties of the external fields which is closely related to the entanglement enhancement. However, there are many 
problems to be solved on the general relationship between the external fields and the entanglement enhancement. 
In particular, we have not shown the properties of the multipartite entanglement. It is obvious that the external 
fields also have effect on qualitative behaviors of the multipartite entanglement; the most famous one is the quantum 
phase transition of the transverse Ising model [24]. In future, we plan to investigate the general properties of the 
enhancement of the multipartite entanglement by external fields. 
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Appendix A: In the case hiop — hivop = 0{f3 with k < in the case (b) 

Here, we discuss the case of /iiop — /ijvop = 0(/3~") with k < in the case (b). In this case, we cannot consider 
the unperturbed states Itiiw) ^ iV'mcdia) liitw) <8) IV'media) independently because their eigenvalues are almost 
degenerate. Then, the magnitudes of the elements Fi and F2 can be different from the ones in (15). We can still 
apply the same calculation to the other parameters {P^J P^^ as in the case k > 0; they are of order ^2ki+2kjv^ 

We here prove that and are of order ^2ki+2kjv qj. higher. In order to prove this, we separate Hi + Hm as 
follows; 

Hi + Hjq = i^LO + SHi^o, 

Hlo = 2(^10? + ^Wop + ho)ai + -(/iiop + /iwop - hQ)a^ 

-{hiop + /iWop)(|titw>(titw| - Ililw)(lil7v|) - /io(|tiijv)(tiijv| - |iltN)(iltjv|), 
^^^LO = —{hiop — h^op — /io)(|tiiw)(tiiw| — lilt n) {lit n\), 

= ^(/llop - ^JVop - ^o)(cti - 0-^), (Al) 

where we define as — 0{P'^°) with < < 1, and regard j35Hi^Q as perturbation. The unperturbed density 
matrix is given by 



P^^l = e-^^''"S (A2) 



where -fftot is defined by 



-^tot — -f^couplc + ^^mcdia + ^LO , (A3) 

where i?mcdia and Hcoupic are given in (19). Because /iq — 0{P^'^°) with < kq < 1, the magnitudes of the 
unperturbed elements {F^\F!f ^} of p^^^ are given by 

^ and 0(/3''i+««-Hko-)^ (^4) 

where 

kq — min(Kjv, 1) and Kq = min(KAr, kq, !)■ (A5) 
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The density matrix in the first-order perturbation is given by 



Ptot 



rtition fund 

and F2 are given by 



1 ' 



+ /? e-^-^- (5iJLoe-^(^-^)*- dx) , (A6) 



where Z^^"* is the partition function of e while + (5Z is the partition function of e The elements Fi 

-Fl = triAr(tl|Ar|ptot||ltw), 

-F2 = triTv (titw IPtot llilw) , (A7) 

where tri^r denotes the trace operation on the spins except the focused spins 1 and N . The first-order perturbations 
of the elements Fi and F2 are 

-(tiljvle-^"^-' |;itw)aitiv|e-'^(i-^)^- llitw)^ (A8) 

and 

-/3(/iiop - /J^op - /io)triAr ^ dx(^(titjv|e-^"^-*' Itilw) (tilivle-^^'-"^^-' IIUjv) 

-(titjv|e-^"^'-|;itAr)aitA^|e-'^(i-^)^'-|;i;jv>), (A9) 

respectively. 

In order to estimate the order of the first-order perturbations of Fi and F2, we introduce {|'0tot'^''^)} a-nd {E^'^l'^^'^} 
as defined in Eqs. (23) and (24), but for H^^l. First, we have 

2«-2 



ri=l 



From the calculation in Ref . [201 , we obtain 



^l^^^t^i^^^^i,] =0(/3--^-). (All) 



Therefore, we have 



where 



Similarly, we obtain 



^(^^||(tii^|e-'^(^-^)^-|;it^)|hO(r^+^"), (A12) 



Z(o)(x)EEtr(e-'3^^-). (A13) 



^(^^ 1 1 (titjv le-''^^- Itiliv) 1 1 = 0{P^- ) , 
^(o)J_^.) 1 1 (tilivle-'^'^-^)^- \UiN) 1 1 = ©(Z?"^ ), 

1 1 (titjv le-''^^- liitiv) 1 1 = 0(/3''^ ), 
^(5j^^3^||(;itA^|e-'^(^-^)^-|;itAr)|| = 0(r"). (A14) 
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As a result, we obtain the first-order perturbations of Fi and F2 as 

0(^«i+K„+i-Ko)^ (A15) 

where we utihzed /3(/iiop — fiNop — ho) = 0(/3^^"") and 

Wc can similarly calculate higher-order perturbations of Fi and F2 to see that they are of order higher than (A15). 
Thus, in the case hiop — /ijvop = 0{/3^'^) with k < 0, F^ and F| are of order higher than {P^J P^J ^^}. 
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Supplementary material for the article "General theory of entanglement enhancement 

by external fields in spin chains" 

In the present supplementary material, we give detailed calculations which cannot be provided in the main text of 
the article "General theory of entanglement enhancement by external fields in spin chains" . 



S.l: Calculation in the case (b) 



Here, we derive the approximated form of j-^JJ,"'^^) in Eq. (23) in the main text, to obtain Eqs. (36)-(39). In order 
to calculate perturbation of Ititjv) IV'media)' employ the general perturbation theory, 

IV') = H'^'^) + (/tot - l^^°^)(V'^°^l) ^(o)^g(o) '^g|V^)' 

where H^^'f is the unperturbed Hamiltonian, Ii/j^")) and iJ^^-* are the unperturbed eigenstate and the unperturbed 

eigenvalue, respectively, and SH is the pcrturbativc Hamiltonian. We calculate each element of j-^^JI.'^^), namely 

{iV'modil,?)} ™ Eq. (23) with ( =tt,tl,lt,U- Wc first calculate iV'mcdiL.it'*' "^^^"^^ ^ivcn by the first-order pertur- 
bation; 

l&Va4t) = l^ltivXIltivl ® /media 



X / 



[/tot - (Itltiv) ® |&ia» ((tltA^I ® (Cedial)] j;(0) ^ g(0) -^couploltlt^) ® l&ia) 



-)N -2 



Y-||.x„|,„ , {{MnI (V'n.cdial)gcouplc(ltltiv) ^ ICedia)) .^ 
„=1 -^modia - "-lop - n-Nop - l-timcdia + "-lop - n.Nop) 

where we put to E^"^"^^.^- hi^p- hNop, H^°^ to //i + //^ + //„,edia, to Ititiv) ® ICcdia) and 6H to //eoupic m 

Eq. (S.l). As has been stated, we assume Ziiop = 0{/3^'^'-) and h^op = 0{j3^^'^) with ki, kjv > 0. On the other hand, 
the eigenvalues {/^niediai order /3° because the media-spin Hamiltonian //media is fixed. Therefore, Eq. (S.2) 

can be approximated by 



2 



media, 4.t' 



E lilt.) ^ IKedia) ^ <-^^-<^'-')^7-('^^^-> ^ 1-^:^"-^^-)) (l + own) . (S.3) 

We can sum the leading term over the label n to obtain 

IV'medil.^t) = r^(liltiv)ait7v| /media) (//eouple |tltiv) ® |&ia) + 0(/3''^)) , 



(S.4) 



where 

2iV-2 

^ ^ IV'modia) ("^medial = /media (S.5) 
ri=l 

is the identity operator in the whole space of the media spins. Using Eqs. (21) and (34) in Eq. (S.3), we have 

iCediL^t) = ® (75-0 ii2) i^i;!) itiv-i) + 7ino 112) iv;ti;>i;A^-i> 

+ "«olt2)i^iI)ltw-i) + ^«„olt2)i^ii)i;w-i) + o(r^)). (S.6) 
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Second, the leading term of IV'mo'dia Ti^ similarly given as follows: 

I V^media.ti' 

N-2 



n=l 
1 

— 2/lArop 
J 



'V- It I \^U,n ^ ((tl^^-l ® (V'mcdial)gcouplc(ltlt^) ® ICnedia)) f, , ^( 
> Itliw) «) IV-modia) ^1 1 + C>{P 

^ -zn Nop ^ 

^mcdia 



-2/i 



Wop 



- 7"no Ii2) IV^il) \iN-l) + Wn. N -l) + 0{P'^-)) . (S.7) 



Third, we calculate IV'me'dia u)' '^hich is given by the second-order pertm-bation; 



n.n' — l 

{{iliN\ » (V'J^ediaD^fcouplcdtUjv) ® K'edia)) ((tlljvl ® «n'odia I) ^^couplc ( Itltjv) ® ICcdia)) 



-^media " '^lop " ^Nop - (^^mcdia + ^lop + ^Wop) ^Todia " ^lop " ^^op " (-E,"^cdia ^ '^lop + ^Wop) 
, ((iliwl ® «,edial)gcouple(|iltiv) ^ IV'ldia)) (UltA^I ^ (V'g.'cdial ) ^couple ( Itltiv) ^ IC°edia)) 



-^modia ^ '^loP ^ ^Nop - {E^^cdia + ^^op + h-Nop) £^,"odia ^ ^^lop - ^Wop " (^"lodia + ^^lop " ^Wop) 



(S.8) 



By utilizing the assumptions hiop = 0{l3~'^'-) > 0, ft,jvop = Ol/?""") > with ki > kat > and i?J^odia = 0(/3''), we 
can approximate (S.8) as 



media/ 



(gilATl (V'^^.ediaDgcouplcdtlliv) ^ ICedia)) 
— 2(ft,iop + tlNop) 





— 2/ijvop 




{{IiIn\ c 


S «cdial)^couplc(|iltiv) ( 


^ l^mcdia)) 




— 2(/liop + /iJVop) 




(Uitivl c 


S (^mcdial)^fcouplc(ltltiv) ( 


» l&ia)) 


— 2/iiop 



(S.9) 



We can calculate the leading terms of (S.9) as 

(Ill|jv)(ll|jv| <^ /mcdia)gcouplc(|tllw)(tliw| /media) ^couple (Itltw) <^ K°cdia)) 

4/iArop('iiop + h^op) 

_^ (|iliiv)(il4-jv| /media)gcouplc(|iltAf)(iltAf| /media) ^couple ( Itltjv) ICcdia)) (-g 

4/liop(/llop + hNop) 

Upon using Eqs. (21) and (34), the numerator of the first term in (S.IO) reduces to 

(|lllw)(lllw| ® /mcdia)i/couplo(|tllw)(tUw| ® Anodia) i/couple ( Itltw) ® ICedia)) 
= (|iliw)(iliAr| ® /media) //couple^ 

X Itllw) ® ilSnJh)\i'll)\iN-l)+tnon2)Hlt)nN-l)+lU^o^^^^ 

^J'IUIn) ® {l's„,\U)\4'll)\iN-l)+ltao\i2)[4^^^^^^^ (S.ll) 
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Similarly, the numerator of the second term in (S.IO) reduces to 

(llllwXIllwl ® /mcdia)i?couplc(|;itiV>Ultw| ® ^mcdia) i^couplc ( |tl Tat) ® K°cdia)) 
= {\UiN){iliN\ ® ^mcdia } ^couple J 

X \Mn) ® {lSno\i2)Hll)nN-l)+lt.,„\i2)Hii)\iN-^^ 

= j'\UiN) ® il's,JUMl)\iN-l)+ltn,^i2)\^^^^^^^^ (S.12) 

which is equal to (S.ll). As a result, we arrive at 

ICedVa.u) = 4/,,^pl^^p (^'^"0 112) I V'll) + 7^»o 112) I^Iji) It^^l) 

+ 7"n„|t2>|^iI>|;iV-l> +^«„„|t2>|^ii>|tiV-l> +0(/3'^")). (S.13) 

Next, we calculate each element of \ip"^[^^), namely {IV'tot c^)} ^q. (23) with ^ =tt,t|,|t,||. By similar calcu- 
lations, we obtain 



ICe'lii.tt) =;jT^(^nolt2>|V':i>ll^-l) +7i«olt2>|V':i>ltAr-l) 



J 

"n„i;2)i^il)i;w-i) +7^^'noi;2)i^ii;)itiv-i) + o(r" 



i&Ya,u) =r^(^^"oi;2)iitI)itAr-i) +7in„i;2)k'i;i)i;A.-i> 

+ «„„it2)iV'iI)itiv-i) + ^i'„„it2)iV'it)i;iv-i) + o(r^)). (S.14) 



The element IV'me'dia 4,t) ^® calculated from 
2W-2 

ICeJii,it)= E lili^) ® ICcdia) 
7i,n' — l 

" (UltTvl ® (CediaD^couplcdtltiv) ® l^^^edia)) ((tltAr| ® (V'ldial ) ^^couplc ( |tli^> ® KCdia)) 



^modia - ^loP + ^^OP - (-^mcdia + ^lop " ^Nop) ^^modia " ^lop + h-Nop ~ (-E^mcdia ~ ^lop ^ ^JVop) 
, (UltA^I ® (V'^J.cdiaD^couplcdllliv) ® |V'r„'cdia» (Uliwl ® (V^tdia I) ^^couplc ( | tliiv) ® ICcdia)) 



^mcdia ~ ^loP + ^^OP ~ (-^mcdia + '^lop " ^Wop) -E^midia ^ ^lop + ^^Vop - (^^media + ^^lop + ^Nop) 



(S.15) 



By utilizing the assumptions /iiop = 0(fi '^^), /lAfop = 0(/3 with k\ > kn > 0, /iiop — ft-JVop = 0(/3 with k > 
and -E'l'Jiodia — we can approximate (S.15) as 



media/ 



(Ult^l ® (V'^^.ediaD^couplcdtltiv) ® IV-^dia)) 



-2(/liop — flNop) 





— 2/lArop 






^^iV-^^ediaD^^couplcdiliivlS 


5 IC'edia)) 




— 2(/liop — tlNop) 




(Uilwlc 


(V'ldial)^couple(ltlijv) S 


^ I'^^media)) 


-2/iiop 



o(r) 



o(r 



(S.16) 



We finally arrive at 

ICediit) -jrA, (7.'^nol;2)|V~:i)liiV-l) +7't„„|;2)K\tt)lt^-l) 

^ftiopftWop ^ 

+ Uno lU) \iN-l) + IWn, \U IV^i!) It^-1> + 0(/3'^" ) + 0(/3«)) , (S.17) 

where k' is defined in (16). 
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S.2: Calculation in the case (c) 

Here, we derive the perturbed form of |0"ot'^) i^- Eq. (63) to obtain Eq. (66). Using Eq. (S.l), the leading term of 
l^totJO Siven by 

Uno,U 'v-'n N^un , (Ull ^ (CediaDgcouplcdtl) ® |Cedia)) , . 

motj - l^l) ® I'Amcdia) p'„o _ ^Tf^^ Tl ^ ' ^^-^^ 

n— 1 media 

"lop l-^mcdia "lopj 

where /iiop = 0{13~'^^) with ki > 1, but {£',j"cdia} ^^'^ of order Then we obtain 

Uno,t\ II \^lAn (Ull <^ ('/'media I) -^couple (Itl) '8' I'/'media)) /^-| , o^ 
I'/'tot.i) = 2^ 111) ® l-Pmcdia) [^ + 0{(3 1)J. (S.2) 

We thereby obtain |(?!>^ot'^) 

l-^t""*^) = Itl) ® l&ia) + ^5^lil> ^5 (l^L\U)\&tN)+<JUMt)\lN) 

— Zftiop ^ 

+ < \t2m:j\tN) + < it2) i0it) 1^^) + oifs""' )) ■ (s-3) 



S.3: Proof for the existence of the maximized entanglement in three-spin systems 

In this section, we prove the statement in Sec. IV. A that the maximized entanglement always exists in the systems 
with the Hamiltonian (81) in the high-temperature limit /3 — >■ 0, except the case of the Ising chain with /i^edia ^ 
We prove it by showing that the entanglement exists by letting and 

hl^hl^ /io(/3), (S.l) 

where we assume 

/io(/3) = 0(/3"''°) with Ko < -1- (S.2) 

We regard the interaction Hamiltonian as perturbation and calculate the leading order of the elements 
{i^i, i^2,Pttj-Pti'Pit'P4-4-} density matrix. Then we show the inequality (8), which is a necessary and sufh- 

cient condition for the existence of the entanglement. 
First, we separate the Hamiltonian (81) as follows: 

Htot = Hi,o + iJint, (S-3) 

where 

Hint = ( J"^T?^T2 + ^''^l^^D + iJ"<^M + J"^l^l) 

= {^J[<^1<^2 + <^l<^2 + li'^l'^2 + '^l'^2)\ + -^[c^^CTa" + <^2'^t + li'^t'^t + <^2'^5)\ (S-4) 

^mcdia arbitrary value, = J(l + 7)/2, = J(l — 7)/2, and we consider ifint as perturbation. Because 

/3/io(/?) — > oo as /3 — > 0, we only have to consider the ground state and the first excited state as the unperturbed 
states, which are given by 

|tit2t3) and Itilats) (S.5) 

with the corresponding eigenvalues 

ei = -2/io - /imcdia + 2J" and £2 = -2/io + /i^edia " 2 J^ (S.6) 

where we assume /i^g^ia > 0, but the following discussion is also applicable to the case of /if„cdia ^ ^- Other excited 
states have the eigenvalues ei + 0{P~'^°) and hence their thermal mixing can be ignored in the limit /3 — >■ 0. 
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We consider the states (S.5) up to the second-order perturbation of Hint- 



'''' -IhlM + TTTT „„v„/ ^I-Ilt2i3) ) (S.7) 



and 



Itli2t3) - 1^ ^T^ , 9 7^ Iilt2t3 



-Itlt2;3> + , ^r.U./'^.,z ^T77Tlili2i3> ) , (S.8) 



2/^0 - 2/iS,,,i, + 2J~' ' (2/10 + 2J-)(2/^o - 2/if„,,,, + 2J^) 

where ki and fc2 £^re the respective normahzation factors. According to the above expressions, ki and ^2 are of order of 
1 + 0{f3'^'^°). By mixing these two states with the Boltzmann weights e~'^'^ and e"'''^^, we obtain the matrix elements 

Let us first consider the case 7 7^ 1. In this case, the leading terms of Fi and ,Jp\fP\4. are given as follows, 
respectively: 



which yields 



Fi - vmm = + o(/3i+2'^°) > (s.io) 



in the limit /3 — > 0. This gives a non-zero value of the concurrence (7) of order /J^"" in the case 7 7^ 0. 

Next, we consider the case 7=1. Because we assumed > > , the equality = = is satisfied in the 
case 7 = 1. In this case, we have to take higher-order approximation because the first term of Eq. (S.IO) vanishes. 
The expansions of F2 and -y/Pti-Pit given, respectively, by 

^'-^0 8 ^ 4 14+^^^ 

Vmm - ^ + ^^"-y^''^' ^ - ^'"^-'f -^ + o(/3^+^-), (s.ii) 



which is followed by 



F2 - = ^^'(/imedia)' + Oi/3^+'-"). (S.12) 



If the media field h^^^^-^i^ is equal to zero, the entanglement vanishes for any values of the local fields hf and h^; this 
is consistent with Theorem 1. For h^^^^^ ^ 0, the concurrence (7) is of order ^i+s^^o jj^ ^j^g ^g^gg ^ = 1 and increases 
as the media field /i^cdia increased. We have thereby proved that the entanglement in a three-spin system with 
Hamiltonian (81) has always a non-zero value except the case of the Ising chain with ft-fnedia = if we choose the local 
fields properly. 



S.4: Calculation of the critical temperature 



Here, we calculate the critical temperature analytically for the use in Sec. V. C. In order to perform the calculation, 
we utilize the numerically confirmed fact that the maximizing local fields take the form /iiop = /i4op ^ 00 in some 
region of the 'Symmetric' phases in Fig. 3. Let us now derive the critical temperature in the limit ft-iop = /i4op ~^ 00. 
In order to calculate the entanglement in this limit, we regard the coupling Hamiltonian i?coupic between the focused 
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spins and the media spins as perturbation. Because of the condition hiop = /i4op oo, we consider only the following 
four unperturbed states: 

Itit4> ® IV^i), |tit4> ® IV'2), |tit4) ® 1^3) and |tit4> ® 1^4), (S.l) 

where {\ipn)}n=i are the eigenstates of the media spins. The states {\ipn)}n=i and their corresponding eigenvalues 
are given by the solution of the bipartite XY spin chain: 

iv-i) = ^(It2;3> - i;2t3», ei = -r - jy, 
= ^(It2;3> + i;2t3», £2 = J" + J^ 

|^3)=5-|;3;4>-5+|t2t3), e3 = -^A{K^^^,J^ + {J--Jy)^ 

1^-4) = 5+|;3;4> + 5_|t2t3), £4 = v/4(/lLdia)' + " •^^)', (8-2) 

where we define 5+ and 5_ as 



1 / 9h^ 
S± = ^Jl± , (S.3) 

Next, we calculate the perturbation of {|tit4) ® |^n)}t=i as in the proof of Theorem 2. The perturbation contri- 
bution of |tit4) ® \ipno) to the leading term of the element Fi is given by 

4ftlop"4op ^ ^ 

where {s„,t„ , li^j^ , Wfi } are defined as the coefficients of {\ipn)}t=i in (S.2): 

IV'n) = S«lt2t3) + tn|t2;3> + Wn|;2t3> + U'n|;2;3>- (S.5) 

Equations (S.4) and (S.5) correspond to Eqs. (40) and (21), respectively. Similarly, the perturbation contribution of 
Itit4) iV'no) to the leading term of the element F2 is given by 

J2 

{I'^WnoSno + T^no^rio + itngUno + S„(,Wno), (S.6) 



4/liop/l4 



op 



which corresponds to Eq. (41). On the other hand, Eqs. (37)"(39) give the perturbation contribution of |tit4) iV'no) 
to the leading terms of the elements {ptJ-iPj-t'^^ii 



/ 2 2 ,4-2 I 22 t 2 \ 



4/i4op 



p 

J4 



respectively. 

We have now all the ingredients for calculating the matrix elements {Fi, F2,p-\^,p^i,Pi^,Pii} in the limit hiop = 
hiop 00. We consider the mixture ofthe perturbed states of |tit4)®l'0i)i Itit4)® l'02), |tit4)<X'l''/'3) and |tit4)'8i|V'4)- 
These four states mix with the Boltzmann weights of e-''(^i+''^i), e-''('2+5e2)^ ^-Pies+Se^) ^^^^ g-/3(e4+<5«4)^ where we 
define the energy perturbations as {Sei}f^i. In the limit /iiop = ^4op — > 00, we have to consider only the leading terms 
of {Fi, F2,p-i-f,p^i,p_i-f,pii}. Therefore, we ignore the energy perturbation {Sei}j^i. We thus arrive at the elements 
of the density matrix as 

ZtotFi - ^[-45_5+7sinh/3|e3|- (l + 72)sinh/3|ei|], 

J2 

ZtotF2 = — 2[-25-5+(l+72)sinh/3|e3| -27sinh/3|ei|], (S.8) 
4/Iq 



29 



and 



ZtotPn = ^tot = 2cosh/3|ei| +2cosh/3|e3| 

7 ' 



{j'S'i+Sl)e-^'^ + {S'i+j^Sl}e-f''-' + (7^ + l)cosh/3|ei| 



ZtotPu 



'0 
J4 



(S.9) 



where we defined /iq = /iiop = h^op with /iq ^ 00. 

By utilizing these parameters, we obtain a necessary and sufficient condition (8) for the existence of the entanglement 
in the form of the following two inequalities: 



Zl,{Fi ~ p^^pu) = [45„5+7sinh/3|e3| + (1 + 7') sinh/?|ei|]' 

- Ztot [h^S^ + Sl)e-^'' + {Si + j^SDe-f^'^ + 2^ cosh/?|ei|] > 0, 



(S.IO) 



or 



Zlot{Fl-VnVi^) = [25_5+(l +7')sinh/3|e3| +27sinh/3|ei|] 
- [(7252 + S\)e-^'' + {Si + 7252 )e"^^^ + (72 + 1) cosh /3|ei|] > 0. 



(S.ll) 



